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Abstract 

This  project  is  focused  on  reduced  basis  approximation  methods,  associated  rigorous 
and  sharp  a  posteriori  error  bounds,  and  offline-online  computational  strategies  for  the 
rapid  and  reliable  solution  of  parametrized  elliptic,  parabolic,  and  more  recently  hy¬ 
perbolic  partial  differential  equations  relevant  to  mechanics  from  the  quantum  through 
the  meso-scale  to  the  macro-scale.  Typical  equations  and  applications  of  interest  include 
Density  Functional  Theory  for  solid  state  property  calculations,  the  Boltzmann  equation 
for  microscale  gas  flows,  the  Navier-Stokes  equations  for  natural  convection  calculations, 
elasticity  for  stress  intensity  factors/brittle  failure,  and  Helmholtz  and  the  wave  equa¬ 
tion  for  acoustic  waveguide  applications.  Of  particular  interest  is  real-time  and  robust 
parameter  estimation  with  application  to  detection,  nondestructive  evaluation,  adaptive 
design/optimization,  and  control. 

In  the  online/deployed  stage,  we  can  provide  results  for  key  engineering  outputs  in 
real-time  without  loss  of  accuracy  or  reliability:  the  outputs  provided  —  in  milliseconds 
(online)  —  by  our  approach  are  provably  indistinguishable  from  the  outputs  provided 
—  typically  in  many  minutes  or  even  hours  —  by  classical  methods. 

Our  web  site  http :  //august  ine .  mit .  edu/  contains  an  interactive  MATLAB®  demo 
(from  the  home  page  click  on  worked  problems),  as  well  as  a  summary  of  the  methodol¬ 
ogy  (from  the  home  page  click  on  methodology)  and  a  compendium  of  our  publications 
(from  the  methodology  page  click  on  Technical  Papers). 

Introduction 

Engineering  analysis  requires  the  prediction  of  a  (or  more  realistically,  several)  “out¬ 
put  of  interest”  s'Gt  —  related  to  energies  or  forces,  stresses  or  strains,  flowrates  or 
pressure  drops,  temperatures  or  fluxes,  intensities  or  radar  cross-sections  —  as  a  func¬ 
tion  of  an  “input”  parameter  P-vector  p  €  V  C  Rp  —  related  to  geometry,  physical 
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properties,  boundary /initial  conditions,  and  loads/sources. 

These  outputs  se(/x)  are  often  functionals  of  a  field  variable 

se(v)  =£(ue(n)\fi)  ,  (1) 

where  xxe(/x)  £  Xe  —  say  displacement,  velocity,  or  current  —  satisfies  a  /i-parametrized 
partial  differential  equation 

a(ue(/j,),v;ij,)  =  /(u;/x),  V v  £  Xe  .  (2) 

Here  Xe  is  the  appropriate  function  space,  and  a  (respectively  £,  f)  are  continuous 
bilinear  (respectively,  linear)  forms.  (We  consider  for  simplicity  in  this  brief  exposition 
the  linear  elliptic  case.) 

In  general,  we  can  not  find  the  exact  (our  superscript  “e”  above)  solution,  and  hence 
we  replace  se(/x),  ue(/x)  with  a  (say)  Galerkin  finite  element  approximation, 
given  /x  £  T>, 

s^in)  =  i(uN(n))ix)  ,  (3) 

where  £  X ^  satisfies 

a(uM{v),v\n)  =  /(v;/x),  V  v  £  X*  .  (4) 

Here  X*  c  ATe  is  (say)  a  standard  finite  element  approximation  subspace  of  dimension 

a r. 


Unfortunately,  to  achieve  the  desired  accuracy,  M  must  typically  be  chosen  very 
large;  as  a  result,  the  evaluation  /x  — >  is  simply  too  costly  in  the  many-query  and 

real-time  situations  often  of  interest  in  engineering.  Low-order  models  —  we  consider 
here  reduced  basis  approximations  —  are  thus  increasingly  popular  in  the  contexts  of 
engineering  analysis,  parameter  estimation,  design  optimization,  and  control. 

In  the  reduced  basis  approach  [1,  5,  14,  16],  we  approximate  sr(ij),vr(n)  —  for 
some  fixed  sufficiently  large  “truth”  M  =  A/t  —  with  sw(/x),  «jv(/x):  given  /x  £  T>, 

sn(h)  =  £{u m(/x);/x)  ,  (5) 

where  u^{n)  £  satisfies 

a{uN{n),v\n)  =  f{y\n),  Vv£WN.  (6) 

(For  simplicity  here,  we  consider  a  purely  primal  approach;  in  actual  practice,  we  pursue 
a  primal-dual  formulation.) 

Here  is  a  problem-specific  space  of  dimension  N  <C  A/t  which  focuses  on  the 
(typically  very  smooth)  parametrically  induced  manifold  of  interest:  \ji  £  V}. 

For  example,  in  the  Lagrange  reduced  basis  approach  (that  we  pursue),  is  the 
span  of  snapshots  1  <  n  <  N,  for  selected  /xn  £  V ;  alternatives  include  the 
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Taylor,  Hermite,  and  “POD”  spaces.  Properly  selected  —  as  we  discuss  below  —  these 
spaces  provide  very  rapid  convergence  u/v(/i)  — » ►  and  hence  Sn(p)  — >  p )  as  N 

increases  [5,  12]. 

The  dramatic  dimension  reduction  afforded  by  these  reduced  basis  spaces  can  be 
transformed  into  very  significant  computational  economies:  the  crucial  additional  in¬ 
gredient  is  an  offline- online  computational  procedure  [2,  10,  11,  17]  that  exploits  affine 
or  approximate  affine  parametric  structure  to  greatly  reduce  the  marginal  cost  of  each 
additional  output  evaluation.  This  strategy  yields  very  large  savings  in  the  many-query 
and  real-time  contexts:  online  (deployed)  complexity  depends  only  on  N  —  typically 
orders  of  magnitude  smaller  than  Aft  —  and  not  on  Aft-  (Of  course,  the  offline  cost 
remains  significant:  some  of  the  Accomplishments  of  this  effort  address  control  of  the 
offline  computational  burden.) 

Our  focus  is  the  development  of  a  posteriori  error  estimators  for  reduced  basis  ap¬ 
proximations  [11,  13,  17,  20]:  inexpensive —  marginal/ online  complexity  independent  of 
Af  —  and  sharp  error  bounds  Ajy(p)  that  rigorously  bound  the  error  in  the  output: 

(7) 

Absent  such  rigorous  error  bounds  we  can  not  determine  if  N  is  too  small  —  and  our 
reduced  basis  approximation  unacceptably  inaccurate  and  decision-making  compromised 
—  or  if  A  is  too  large  —  and  our  reduced  basis  approximation  unnecessarily  expensive 
and  real-time  response  compromised.  In  more  pragmatic  terms,  without  error  bounds 
we  can  not  rapidly  and  rigorously  determine  if  critical  design  or  safety  conditions  and 
constraints  are  satisfied:  for  example,  if  Sn(/j)  represents  a  stress  intensity  factor,  does 
approximate  feasibility /safety  s^{p)  <  C  imply  “true”  feasibility /safety  (/z)  <  (7? 
(In  the  nonlinear  context,  error  bounds  are  crucial  in  establishing  the  very  existence  of 
a  “truth”  solution  u^1  (p)  [9,  19].) 

In  fact,  absent  rigorous  error  bounds,  we  can  not  even  construct  an  efficient  and  well- 
conditioned  reduced  basis  approximation  space  W^.  In  particular  for  higher  parameter 
dimensions  P,  error  bounds  play  a  crucial  role  in  “greedy”  sampling  procedures  [7,  8, 
13,  20]:  we  explore  a  very  large  subset  of  the  parameter  space  V  to  successively  find 
the  points  with  the  largest  (inexpensively  computed)  error  bound ;  we  only  evaluate  the 
truth  solution  for  these  winning  candidates  —  the  “snapshots”  vr%  (p„)  on  the  parametric 
manifold  the  span  of  which  defines  WN.  Conventional  sampling  techniques,  which  do 
not  exploit  error  bounds,  either  can  explore  only  very  small  subsets  of  V  or  require 
prohibitively  expensive  offline  effort.  Effective  sampling  procedures  are  essential  to  rapid 
convergence  and  efficient  representation. 

Accomplishments 

We  describe  here  some  specific  accomplishments, 
a.  Number  of  parameters . 

We  have  made  considerable  progress  in  at  least  two  directions. 
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(i)  Inf-Sup  Lower  Bound.  We  continue  to  improve  the  inf-sup  lower  bound  procedure 
crucial  to  our  a  posteriori  output  bounds  for  difficult  (in  particular,  noncoercive) 
problems.  We  can  now  reasonably  treat  (say,  Helmholtz)  problems  with  several 
resonances  and  0(5)  parameters.  We  have  developed  two  approaches. 

The  first  approach  is  a  higher-order  expansion  based  on  a  new  “natural  norm.” 
The  technique  has  now  been  applied  in  the  coercive,  non-coercive,  and  also  non¬ 
linear  contexts.  (See  S  Sen,  K  Veroy,  DBP  Huynh,  S  Deparis,  NC  Nguyen,  and 
AT  Patera,  “Natural  norm”  a  posteriori  error  estimators  for  reduced  basis  ap¬ 
proximations.  Journal  of  Computational  Physics,  217  (1):  37-62,  2006.  doi: 
j  .  j  cp .  2006 . 02 . 012.) 

The  second  approach  is  a  “Successive  Constraint  Method”  which  converts  the  para¬ 
metric  inf-sup  lower  bound,  via  a  Rayleigh  quotient,  to  a  Linear  Program.  The 
method  offers  good  theoretical  structure,  simplicity,  and  (in  tests  to  date)  quite 
good  offline  and  online  performance.  Future  work  is  required  to  develop  an  a  priori 
convergence  theory,  more  efficient  constraint  selection,  and  faster  Linear  Program 
solution.  A  paper  on  this  subject  has  been  submitted  to  CRAS  Mathematics  (see 
http://augustine.mit . edu/methodology/papers/aptCRAS2006preprint .pdf 
for  the  preprint). 

(m)  Greedy  Sample  Selection.  On  the  theoretical  front  (work  with  A.  Buffa,  Y.  Ma- 
day,  C.  Prud’homme,  G.  Turinici),  related  to  approximation,  we  have  shown  that 
if  a  sequence  of  approximation  spaces  exist  that  provide  sufficiently  rapid  (expo¬ 
nential)  convergence,  then  our  greedy  sample  selection  procedure  will  efficiently 
find  an  “associated”  space  that  also  provides  (somewhat  less)  rapid  exponential 
convergence.  Hence  although  we  do  not  yet  know  a  priori  when  to  expect  rapid 
convergence,  we  do  know  (roughly)  that  if  in  theory  sufficiently  fast  exponential 
convergence  is  possible,  then  in  practice  our  procedure  will  also  exhibit  very  good 
approximation  properties.  A  paper  is  in  preparation. 


b.  Generality  of  parametric  dependence. 

We  have  earlier  developed  a  technique,  the  “Empirical  Interpolation  Method”  [3], 
which  permits  us  to  develop  efficient  reduced  basis  approximations  for  problems  with 
non-affine  parameter  dependence  [6]. 

In  this  project  we  have  developed  completely  rigorous  error  estimators  for  non-affine 
problems:  these  estimators  rely  on  a  high-order  Taylor  series  development  and  on  the 
demonstrably  rapid  convergence  of  reduced  basis  approximations  for  parametric  deriva¬ 
tives  of  the  field  variable.  We  have  demonstrated  the  error  bounds  for  several  simple 
model  problems.  Future  work  will  consider  more  ambitious  problems  with  non-affine 
geometric  variation. 

c.  The  Helmholtz  equation. 

We  have  developed  reduced  basis  approximations  and  a  posteriori  error  estimators  for 
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Figure  1:  Acoustic  Low-Pass  Filter  Geometry. 


waveguides  with  accurate  outflow/radiation  conditions  for  multiple  propagating  modes 
(and  of  course  also  evanescent  modes).  These  methods  will  permit  us  to  consider  a  large 
range  of  acoustic  waveguide  problems  related  to  noise  mitigation  or  detection  (say,  of 
buried  objects  in  the  seabed). 

As  an  example,  we  consider  here  the  acoustic  low-pass  filter  waveguide  element  shown 
in  Figure  1.  The  parameters  are  the  (non-dimensional)  frequency  squared,  /xj  =  k2,  and 
the  width  of  the  low-pass  expansion  (relative  to  the  nominal  waveguide  width),  /x2  =  H] 
we  consider  the  parameter  domain  /x  £  V  =  [0.1, 5.0]  x  [1.75, 2.25]. 

The  governing  equation  is  the  Helmholtz  acoustic  equation  with  inhomogeneous  Neu¬ 
mann  (imposed  velocity)  at  inflow,  Fin,  damping  on  the  liner  surface,  Fijner,  waveguide 
radiation/propagation  conditions  (a  Robin  condition)  on  rout,  and  homogeneous  Neu¬ 
mann  (zero  velocity)  on  all  other  boundaries. 

The  output  of  interest  is  a  transmission  coefficient  TC,  the  log  of  the  ratio  of  the 
ouput  pressure  intensity  to  the  input  pressure  intensity. 

We  present  in  Figure  2  the  reduced  basis  output  TCjv  for  N  =  35  as  a  function 
of  /xi  =  k2  at  several  values  of  /x2  =  H.  We  also  indicate  the  error  bars:  the  truth 
FEM  solution  will  with  certainty  reside  within  the  indicated  (rather  tight)  error  bars. 
For  this  particular  problem  the  online  computational  savings  —  the  ratio  of  the  time 
to  calculate  the  truth  output,  ^‘(/x)  to  the  online  time  to  compute  the  reduced  basis 
output  and  error  bound,  Sn{p)  and  A sN(n)  —  is  only  0(5)  rather  than  the  usual  0(100) 
(see  e.  below  and  Table  1):  first,  the  truth  model  for  this  simple  test  case  is  quite  coarse; 
second,  the  Linear  Program  performance  for  our  new  SCM  inf-sup  lower  bound  is  quite 
poor.  The  latter  will  be  the  subject  of  future  work. 

We  are  also  collaborating  with  the  group  of  J.  Hesthaven  of  Brown  University  and 
Y.  Maday  of  University  of  Paris  VI  on  reduced  basis  methods  and  a  posteriori  error 
estimators  for  electromagnetic  Helmholtz  applications. 
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Figure  2:  Low-pass  filter  transmission  coefficient  as  a  function  of  k 2  for  different  H 
(see  Figurel).  As  expected,  for  low  k 2  the  signal  passes;  for  moderate  k2  the  signal  is 
attenuated;  and,  at  higher  k 2  (due  to  harmonics)  the  signal  rebounds. 

d.  Real-time  robust  parameter  estimation  and  uncertainty  quantification  for  detection, 
adaptive  design,  and  control. 

To  date  we  have  further  developed  our  parameter  estimation  techniques,  in  particular 
for  the  efficient  construction  of  “possibility  sets”  through  application  of  various  optimiza¬ 
tion  procedures.  We  have  also  recently  discovered  that,  for  an  interesting  class  of  elliptic 
and  parabolic  problems  (e.g.,  related  to  transient  thermal  detection  of  cracks),  we  can  in 
fact  achieve  complete  rigor  while  still  maintaining  reasonable  computational  efficiency; 
the  critical  enablers  are  the  reduced  basis  output  approximations  and  a  posteriori  error 
estimators.  A  paper  has  been  submitted  to  the  journal  Inverse  Problems. 

e.  Stress  Intensity  Factors. 

Stress  Intensity  Factors  (SIF)  are  crucial  in  the  prediction  of  fatigue-induced  crack 
growth  and  potential  brittle  failure.  In  the  past,  SIF  calculations  are  either  too  crude 
(handbooks)  or  too  expensive  (direct  finite  element  calculation)  to  be  useful  in  the 
most  important  contexts:  many-query  evaluation  as  needed  in  fatigue  studies;  real-time 
evaluation  as  needed  in  embedded/health  monitoring/lifing  analyses. 

We  have  developed  a  new  formulation  for  SIFs  particularly  well  suited  to  fast,  sharp, 
and  rigorous  reduced  basis  prediction  and  a  posteriori  error  estimation.  This  new  for¬ 
mulation,  and  associated  reduced  basis  treatment,  is  described  in  [8].  This  paper  also 


7 


t  t  t  t  t 


1  i  1  1  I 

Figure  3:  Geometry  for  Stress  Intensity  Factor  calculation. 


N 

E  N 

£n 

Online 

Computation  Savings 

5 

1.04E-01 

1.66E+02 

467 

10 

6.01E-02 

8.72E+01 

417 

20 

9.08E-03 

4.39E-01 

248 

30 

2.36E-03 

1.31E-01 

191 

35 

9.42E-05 

3.17E-02 

136 

40 

4.58E-05 

8.86E-03 

125 

Table  1:  Reduced  basis  error,  a  posteriori  error  bound,  and  online  computational  sav¬ 
ings;  the  latter  is  the  time  to  compute  the  truth  solution  divided  by  the  time  to  compute 
(online)  the  reduced  basis  output  and  associated  a  posteriori  error  bound. 

presents  several  examples,  one  of  which  we  briefly  describe. 

We  consider  here  the  geometry  shown  in  Figure  3:  a  pair  of  symmetric  cracks  each 
of  length  d  (measured  from  the  centerline)  emanating  from  a  hole  of  radius  R  in  a  (plain 
strain)  specimen  of  height  2 L;  all  lengths  are  non-dimensionalized  by  the  specimen  width. 
Our  interest  is  in  the  (linear  elastic)  Stress  Intensity  Factor,  SIF,  as  a  function  of  our 
three  parameters:  p\  —  d,  p2  —  R ,  and  ps  =  L.  The  SIF  serves  to  predict  brittle  failure, 
or  fatigue-induced  crack  growth. 

We  show  in  Table  1  the  error  in  the  reduced  basis  Energy  Release  Rate  (ERR  —  the 
square  of  the  SIF),  E^,  the  a  posteriori  error  bound  for  the  reduced  basis  prediction, 
£jv,  and  the  (online)  computational  savings  of  the  reduced  basis  approach  relative  to 
the  FEM  truth;  the  error  results  reported  are  the  maximum  over  a  large  parameter  test 
sample.  We  observe  rapid  convergence  with  N,  relatively  sharp  (and  rigorous)  error 
bounds,  and  very  large  computational  savings. 

The  reduced  basis  results  can  now  be  used  with  confidence  in  predicting  failure  or 
crack  growth  in  real-time,  a  critical  capability  for  component  design,  health  monitoring, 
and  optimal  (and  safe)  mission  planning.  At  present,  no  other  approaches  can  provide 
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accurate,  real-time,  and  reliable  predictions  for  critical  failure  indicators. 

Future  work  will  address  more  complex  loadings  and  geometric  configurator.  If 
successful,  our  work  could  replace  much  of  the  existing  SIF  technology. 

f.  Hyperbolic  problems. 

We  have  made  progress  on  two  classes  of  hyperbolic  problems:  first-order  wave  equa¬ 
tions  —  in  particular  the  Boltzmann  equation;  and  the  second-order  wave  equation  — 
in  particular  relevant  to  time-domain  acoustics  or  elastodynamics  or  Maxwell. 

The  Boltzmann  problem  is  a  very  difficult  equation  that  involves  independent  space 
and  velocity  coordinates  and  (through  the  collision  term)  non-local  interactions;  the 
equation  is  increasingly  relevant  given  the  growth  in  microtechnology.  We  have  devel¬ 
oped  a  Petrov  Galerkin  stabilized  formulation  of  the  Boltzmann  equation  that  in  turn 
permits  us  to  develop  very  rapidly  convergent  reduced  basis  approximations  and  associ¬ 
ated  rigorous  error  bounds  (e.g.,  for  the  flow  rate  through  a  microchannel).  At  present 
we  can  consider  only  one  space  and  velocity  coordinate  and  a  simple  collision  model 
with  two  parameters  (the  Knudsen  number  and  the  accommodation  coefficient).  This 
work  is  in  collaboration  with  E.  Rpnquist;  a  paper  on  this  topic  will  appear  shortly  [15]. 

For  the  linear  second-order  wave  equation  we  have  found,  perhaps  not  surprisingly, 
that  for  sufficiently  smooth  initial  conditions  we  can  develop  rapidly  convergent  reduced 
basis  approximations.  More  surprisingly,  we  have  found  that  we  can  develop  rigorous 
a  posteriori  error  bounds  for  the  error  (in  the  energy  norm)  as  a  function  of  time  and 
parameter.  A  description  of  the  approach  and  some  first  simple  numerical  tests  are 
provided  in  [18]. 

g.  Nonlinear  problems. 

We  have  continued  our  efforts  on  the  Navier-Stokes  equations  and  natural  convection 
(relevant,  for  example,  to  materials  processing  flows);  we  are  also  considering  a  suite  of 
flows  in  canonical  channel  geometries  relevant  to  the  design  of  new  microfluidic  devices. 

On  a  different  front,  Density  Functional  Theory  quantum  treatment  of  crystalline 
(periodic)  solid  materials,  we  have  recently  initiated  our  first  comparisons  with  conven¬ 
tional  (plane  wave/Fourier  spectral)  approaches.  It  appears  that,  for  our  simple  model 
problems,  our  approaches  can  significantly  reduced  the  requisite  degrees  of  freedom. 
This  work,  in  collaboration  with  Y.  Maday,  C.  Le  Bris,  and  E.  Cances,  is  described  in 

[4]- 
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